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The spatial distribution of interacting chemical fields is 
investigated in the non-diffusive limit. The evolution of fluid 
parcels is described by independent dynamical systems driven 
by chaotic advection. The distribution can be filamental or 
smooth depending on the relative strength of the dispersion 
due to chaotic advection and the stability of the chemical 
dynamics. We give the condition for the smooth-filamental 
transition and relate the Holder exponent of the filamental 
structure to the Lyapunov exponents. Theoretical findings 
are illustrated by numerical experiments. 



The Lagrangian description of fluid flows, considered 
in the framework of chaotic dynamical systems, has given 
important insight into mixing and transport M . In many 
situations of practical interest scalar fields describing 
some property of the fluid are not just simply advected 
by the flow, but are active in the sense that they have 
their own dynamics, in general coupled to the transport 
and mixing process. (In the following we will refer to 
this as the 'chemical' dynamics of the system.) Typical 
examples of active fields in mixing flows are concentra- 
tions of reacting chemicals |^,|| (in industrial processes 
or the atmosphere) and interacting biological populations 
of microorganisms in a fluid (e.g. plankton populations 
stirred by oceanic currents 0]). The spatial structure 
of such fields often has complex filamental character . 
Previous work has investigated the temporal evolution 
of reactions like A + B^CorA + B^ 2B in closed 
flows P,PP] as an initial value problem. The absence 
of chemical sources in these cases necessarily implies a 
homogeneous final state of the system. The same reac- 
tions were also studied in open chaotic flows || where a 
stationary fractal distribution arises due to the proper- 
ties of the underlying (chaotic scattering-like) advection 
dynamics M. Here we show that in the case of stable 
chemical dynamics (in a sense defined below) and in the 
presence of chemical sources persistent filamental fractal 
patterns can also arise in closed flows. 

The governing equations for a set of N interacting 
'chemical' fields Cj mixed by a flow v(r,t), independent 
of the chemical dynamics, can be written as 



dd 

dt 



+ v(r, t) • VCi = J=i[Ci(r), .., C N (r),r], (1) 



spatial derivatives of the fields Ci, e.g. a Laplacian term 
representing diffusion. 

Let us assume that the advective transport dominates 
and diffusion is weak. If we neglect non-local processes 
like diffusion on the right hand side of (pi) F{ becomes 
a simple function of the local concentrations and coordi- 
nates. In this case the Lagrangian representation leads 
to a considerable simplification of the equations, leading 
to a low-dimensional dynamical system 



rfr 
dt 



v(r,t), 



da 

dt 



F i (C 1 ,C 2 ,..,C N ,r), %=1,.,N, 



(2) 



where the second set of equations describes the chem- 
ical dynamics inside a fluid parcel that is advected by 
the flow according to the first equation. The coupling 
between the flow and chemical evolution is non-trivial if 
some of the functions Fi depend explicitly on the coordi- 
nate r. In applications this dependence of Fi can appear 
as a consequence of spatially varying sources, spatially 
varying (e.g. temperature dependent) reaction rates or 
reproduction rates of biological species. It is therefore 
natural to include this dependence in the model to be 
considered. 

Although we use a Lagrangian representation, our aim 
is to follow not the chemical evolution along individ- 
ual trajectories, but the spatiotemporal dynamics of the 
chemical fields Ci(r, t) that is equivalent to the evolution 
of the ensemble of fluid parcels under the dynamics (g) . 
Thus the original problem defined in an infinite dimen- 
sional phase space is reduced to an ensemble problem in 
a low-dimensional dynamical system. 

By neglecting diffusion, we may miss certain classes 
of behaviour such as propagating waves, typical for 
reaction-diffusion systems M]. However, as we shall see 
later, we capture a non-trivial behaviour that we be- 
lieve to be characteristic of the full advection-reaction- 
diffusion problem when diffusion is weak. 

We assume that the flow v(r, t) is two-dimensional in- 
compressible and periodic in time with period T . These 
conditions in general lead to chaotic advection even in 
case of simple spatially smooth (non-turbulent) velocity 
fields. Since the advection is independent of the chemi- 
cal dynamics it can be characterised by its own Lyapunov 
exponents Af > Af\ Incompressibility implies that the 
advection is described by a conservative dynamical sys- 



where i = 1, .., N. The operators Ti in general contain tern and thus Af = — \ 2 



For numerical investigations we used a simple time- 
periodic flow that consists of the alternation of two steady 
sinusoidal shear flows in the x and y direction for the 
first and the second half of the period, respectively. The 
flow is defined on the unit square with periodic boundary 
conditions by the velocity field 

v x (x,y,i) = -^0(77 -£modTJ cos(27ry) 

v y (x,y,t) = -—@(t modT- —J cos(27rx) (3) 

where 0(x) is the Heavyside step function. U — > cor- 
responds to an intcgrable limit of the advection problem. 
We will consider the case U — 1.0 producing a flow that 
consists of one connected chaotic region. The param- 
eter T controls the relative time-scale of the flow and 
reactions. By changing T we can vary the Lyapunov 
exponents of the flow without altering the shape of the 
trajectories and the spatial structure of the flow. 

Since the trajectories r(i) are chaotic the chemical dy- 
namics (|2J) corresponds to a chaotically driven dynamical 
system. This subsystem can be characterised by the set of 
chemical Lyapunov exponents Af > ASp > ... > A^, that 
depend on the driving by the chaotic advection. Here 
we will only consider the simplest situation, when the 
largest chemical Lyapunov exponent is negative and, for 
a fixed trajectory r(i), the chemical evolutions converge 
to the same globally attracting chaotic orbit for any ini- 
tial condition in the chemical subspace. In the case of 
no explicit space dependence in the functions Fi (i.e. no 
driving) this would correspond to a globally attracting 
fixed point of the chemical subsystem. 

A simple example of such a system (relevant for at- 
mospheric photochemistry) is the decay of a chemical 
species (N = 1) produced by a non-homogeneous source, 
with chemical dynamics 



C = o(r) - bC, 



(4) 



for which X c = —b. 

First we investigate the temporal evolution of the 
chemical fields. Although the Lagrangian variables Cj(i) 
have a chaotic time-dependence according to the positiv- 
ity of the largest Lyapunov exponent for the full dynami- 
cal system (g) , it can be shown that the Eulerian chemical 
field Ci(r,t) is asymptotically periodic in time with the 
period of the flow. In order to obtain the values of the 
fields at point r at time t one can i) integrate the advec- 
tion backward in time and use the obtained trajectory 
r(i') (0 < t' < t) and the initial values of the chemical 
fields at the point r(0) to it) integrate the chemical dy- 
namics forward in time from to t. The value of the 
field at the same point at time t + nT can be obtained 
similarly. The resulting backward trajectory will be the 
same (due to the periodicity of the flow) but longer. By 
integrating the chemical evolution forward in time for nT 



we obtain a problem equivalent to the previous one with 
a different set of initial concentrations, that according to 
the assumptions made above (Af < 0) converge to the 
same orbit. Consequently the chemical fields are asymp- 
totically periodic in time 

lim Ci(r, (m + n)T + T) = lim C t (r,mT + r), (5) 

where r = t mod T, defining an asymptotic chemical field 
C°°(r,T). 

In the following we investigate the spatial structure of 
the chemical fields. For this we calculate the difference 



SCi = Ci(r + 6r,t) - Ci(r,t), \6r\ < 1, 



(6) 



that can be obtained by integrating (|2|) along two tra- 
jectories ending at the preselected points r and r + <5r at 
time t. The time evolution of the distance |<5r(i')| (for 
t — t 1 ~^> 1) can be estimated from the time reversed ad- 
vection dynamics as 



\Sr{t')\ ~ |<Jr(i)|e 



\ F (t-t') 



(7) 



where A F = Af > for almost all final orientations 
n(t) = Sr(t)/\Sr(t)\. The only exception is the unsta- 
ble contracting direction of the time reversed flow corre- 
sponding to X F — X2 < 0. 

By expanding (||) around the chaotic orbit Cj(r(t)) we 
obtain a set of linear equations 



with initial condition 



^dR 



— -SCj + V r F t ■ Sr 



«5Q(0) = C°(r(0) + 5r(0))-C I (r(0)), 



(8) 



(9) 



where C°(r) are the initial chemical fields. 

In the simplest case N = 1 the solution of (||) can be 
written explicitly as 

SC(t) = V r C° • 5r(0)e xt + f V r F ■ < 5r(i')e A( ^ t ' ) rfi' 

(10) 



where 



1 /"' r)F 



(11) 



that becomes equal to the chemical Lyapunov exponent 
A c in the t — » 00 limit. The first term in (nO) represents a 
deviation due to non-homogeneous initial conditions and 
the second one describes the effect of different histories of 
the two trajectories appearing via the space dependence 
of F. 

Using (fil) we obtain from (FTo) for the projection of the 
gradient of C to a direction n 



SC(t) 

WW\ 



= VC°n(0)e( AF+AC )' + 

+ f VF(r(f'))n((')e (AF+AC)M) (it'. (12) 
Jo 



The convergence of the right hand side of the above equa- 
tion for t — » oo depends on the sign of the exponent 
X F + X c . If X F < \X C \ the convergence of the chemical 
dynamics is stronger than the dispersion of the trajecto- 
ries in the physical space and results in a smooth field 
C°°(r,T). On the contrary, if X F > \X C \ the memory 
of the chemical dynamics decays too slowly to forget the 
different spatial histories (or initial conditions). In this 
case the limit does not exist and the field C°°(r, r) has an 
irregular structure that is almost nowhere differentiablc. 
There exists, however, at each point one special direc- 
tion in which the derivative is finite. This direction is 
the contracting direction in the time-reversed advection 
dynamics corresponding to the negative Lyapunov expo- 
nent of the flow X F . Thus we suggest a precise definition 
of a ftlamental structure, being a non-diffcrentiable field, 
that is still smooth in one direction at each point (with 
that direction itself varying smoothly) . 

The irregular chemical field can be characterised by 
its Holder exponent a defined as 8C(Sr) ~ \Sr\ a , where 
< a < 1 and a = 1 for a differentiablc function. In our 
case if A F + X c > 



SC(t) ~ \6r(t)\e^ F+ * C ». 
Expressing time from (pi) as 



1 |fa(0)| 

A^ ln |<5r(t)l 



and inserting it in (n3|) we obtain 
SC(t) ~ |5r(0)| 



(i+* c /* F )\8r(t)\- x °J xF . 



(13) 



(14) 



(15) 



For very long times and in a closed flow |<5r(0)| will satu- 
rate at a finite value in the backward advection dynamics. 
Thus the Holder exponent is a — \X C \/X F . Certainly, dif- 
fusion would smooth out the small scale filamcntal struc- 
tures below a certain diffusive scale (approaching zero for 
smaller and smaller diffusivities) , setting a cut-off for the 
scaling relation (15). Nevertheless, above the diffusive 
scale filamental structures will persist for arbitrarily long 
time since, in the presence of chemical sources, the effect 
of diffusion is balanced by the continuous generation of 
small scale structures by the chaotic advection. 

As an example we consider the system (Eh, which is 
the simplest possible system that exhibits the smooth- 
filamental transition, using a source term of the form 
a(x,y) = 1 + ai sin (27rx) sin (2ity). Numerically calcu- 
lated chemical fields (obtained by a backward integra- 
tion of the advection problem and forward integration of 
the chemical dynamics pd|) of both types are shown in 



Fig. la and b. Figure 2 shows a section of the smooth and 
filamental fields of Fig.l along the line y = 0.25. We also 
measured the box counting fractal dimension of the func- 
tion C(r) along a cut in the x direction, which is related 
to the Holder exponent by D = 2 — a J12] . Numerically 
computed values agree well with the theoretical predic- 
tion as shown in Fig. 3. If the flow does not consist of 
only one connected chaotic region, but is composed by 
chaotic regions separated by KAM tori, the Lyapunov 
exponents of the flow are different in each chaotic region 
and the structure of the chemical field can be smooth in 
certain regions but filamentary in others (Fig la). 

In the more general case N > 1 the same smooth- 
filamental transition can be obtained. Let us once again 
consider at a particular time t the difference SCi(t) in 
the chemical fields over a small preselected displacement 
<5r. The time dependence of a set of initial deviations 
6Ci(Q) and <5r(0), for the dynamical system (g) has the 
asymptotic form 



N 

X- 

2 = 1 



1/2 



N 



(J2^f(t) + Sr 2 (t)) ~ (£*C?(0) + 5 r 2 (0)) 



1/2 



A/ 



i=l 



(16) 



where A is one of the Lyapunov exponents of the system. 
For a typical choice of the final deviation Sr(t), <5r will 
be divergent in the time reversed dynamics, and conse- 
quently contracting in the forward direction. Thus, the 



contribution of 5r decays in the left hand side of (16) and 
A cannot be the positive Lyapunov exponent X F so the 
typical value of A will be the second largest Lyapunov ex- 
ponent. There are thus two possibilities A = X F = — Xf 
(if X F > Xf ) or A = Af (otherwise). We can now divide 



both sides of @ by Sr(t) using (@) 



N 



Sr 2 (t) 



1/2 



N 

E 



<5r 2 (0) 



1/2 



(Af +A)t 



(17) 



Thus, the chemical fields Cj(r) become non-differcntiable 
in the t — > oo limit if Af > Xf \ . 

Numerically we have also observed the smooth- 
filamental transition in some two-component sys- 
tems (Brusselator autocatalytic reaction model [fl3[ , 
phytoplankton-zooplankton population dynamics model 
|14| ) for parameters that satisfies the condition Af < 0, 
and obtained distributions very similar to that shown in 
Fig.l. 
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FIG. 1. Smooth (a) and filamental (b) distributions of the 
decaying substance (eqMj) mixed by the flow (eq. 0) after 20 
periods using a\ = 0.1 ,b = 1.0 ,U = 1.0 and initial condition 
C(r, t = 0) = 0.0. The integration has been done for 200 x 200 
points with final positions on a rectangular grid. For U = 1.0 
the numerically calculated Lyapunov exponent of the flow is 
A F « 2.35/T. The period of the flow is T = 5.0 (A F < b) in 
(a) and T = 1.0 (A F > b) in (b), respectively, (c) coexistence 
of smooth and filamental structures for U = 0.5 and T = 1.0. 



FIG. 2. Sections of the smooth (a) and filamental (b) fields 
shown in Fig.l along the line y — 0.25. 



FIG. 3. Box counting fractal dimensions of the functions 
C(x) shown in Fig. 2. Number of boxes needed to cover the 
graph of the function C(x) vs. the box size I (squares) and 
slopes corresponding to the relation D = 2 — b/X F . 
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